The gut microbiota-immune-brain axis in a wild vertebrate: dynamic interactions and health impacts

The gut microbiota-immune-brain axis is a feedback network which influences diverse physiological processes and plays a pivotal role in overall health and wellbeing. Although research in humans and laboratory mice has shed light into the associations and mechanisms governing this communication network, evidence of such interactions in wild, especially in young animals, is lacking. We therefore investigated these interactions during early development in a population of common buzzards (Buteo buteo) and their effects on individual condition. In a longitudinal study, we used a multi-marker approach to establish potential links between the bacterial and eukaryotic gut microbiota, a panel of immune assays and feather corticosterone measurements as a proxy for long-term stress. Using Bayesian structural equation modeling, we found no support for feedback between gut microbial diversity and immune or stress parameters. However, we did find strong relationships in the feedback network. Immunity was negatively correlated with corticosterone levels, and microbial diversity was positively associated with nestling body condition. Furthermore, corticosterone levels and eukaryotic microbiota diversity decreased with age while immune activity increased. The absence of conclusive support for the microbiota-immune-brain axis in common buzzard nestlings, coupled with the evidence for stress mediated immunosuppression, suggests a dominating role of stress-dominated maturation of the immune system during early development. Confounding factors inherent to wild systems and developing animals might override associations known from adult laboratory model subjects. The positive association between microbial diversity and body condition indicates the potential health benefits of possessing a diverse and stable microbiota.

The gut microbiota-immune-brain axis is a feedback network which influences diverse physiological processes and plays a pivotal role in overall health and wellbeing.Although research in humans and laboratory mice has shed light into the associations and mechanisms governing this communication network, evidence of such interactions in wild, especially in young animals, is lacking.We therefore investigated these interactions during early development in a population of common buzzards (Buteo buteo) and their e ects on individual condition.In a longitudinal study, we used a multi-marker approach to establish potential links between the bacterial and eukaryotic gut microbiota, a panel of immune assays and feather corticosterone measurements as a proxy for long-term stress.Using Bayesian structural equation modeling, we found no support for feedback between gut microbial diversity and immune or stress parameters.However, we did find strong relationships in the feedback network.Immunity was negatively correlated with corticosterone levels, and microbial diversity was positively associated with nestling body condition.Furthermore, corticosterone levels and eukaryotic microbiota diversity decreased with age while immune activity increased.The absence of conclusive support for the microbiota-immune-brain axis in common buzzard nestlings, coupled with the evidence for stress mediated immunosuppression, suggests a dominating role of stress-dominated maturation of the immune system during early development.Confounding factors inherent to wild systems and developing animals might override associations known from adult laboratory model subjects.The positive association between microbial diversity and body condition indicates the potential health benefits of possessing a diverse and stable microbiota.

Introduction
A growing body of evidence highlights a crucial relationship among gut microbes, the environment and host immunity, which also molds the maturation and function of the central nervous system (CNS) (Fung, 2020;Chin et al., 2020).This intricate feedback system encompasses the brain, the autonomic and enteric nervous systems, the endocrine and immune systems, and the gut microbiome (Martin et al., 2018).Gut-brain communication occurs through various mechanisms, including neuroanatomical pathways and the neuroendocrine axis (hypothalamic-pituitaryadrenal or HPA axis).However, increasing importance has been given to the interactions between the microbiota and the immune system as a fundamental pathway regulating microbiota-gut-brain communication (Fachi et al., 2019).For example, the products resulting from the metabolic activities of the gut microbiota include bioactive peptides such as neurotransmitters, short-chain fatty acids, intestinal hormones, and branched-chain amino acids that play pivotal roles in regulating communication between the gut, brain and the immune system (Sperringer et al., 2017;Cryan et al., 2019;Ding et al., 2020;Krishnamurthy et al., 2023).These metabolites are able to enter the circulatory system to communicate with the brain, eliciting stimulation of the HPA axis (Lynch et al., 2024;Sudo et al., 2004).Additionally, they exert a direct influence on the mucosal immune system, where specific metabolites can function as immune signaling molecules, activating or enhancing systemic immune responses (Ding et al., 2020;Long-Smith et al., 2020).In turn, the brain regulates the gastrointestinal tract and overall organismal homeostasis (Mindus et al., 2021).
The gut microbiota contributes to the development of the host immune system during early life stages, regulates it to maintain gut homeostasis and protects the host from the colonization of potential pathogens (Lo et al., 2021;Berman et al., 2023).As such increasing importance has been given to the neonatal microbiome and its developmental trajectory during early life which affects host immunity and neural activity maturation (Jašarević and Bale, 2019;Ratsika et al., 2023).For example, in humans, the mode of delivery significantly influences the normal colonization of the gut microbiota (Ratsika et al., 2023).C-section deliveries have been linked to immune disturbances (allergies and asthma) (Roduit et al., 2009;Bisgaard et al., 2011) and disruptions in the structure and function of the CNS, with indications of increased neuronal cell death (Castillo-Ruiz et al., 2018) and deficits in early communication skills and social behavior later in life (Morais et al., 2020).In murine models, penicillin treatment during late pregnancy and early postnatal life has been shown to induce alterations in the gut microbiota, leading to heightened cytokine expression, altered integrity of the blood-brain barrier, and significant differences in behavior (Ratsika et al., 2023).Interactions among these functional systems, particularly during sensitive developmental periods, carry consequences later in life, influencing behavior, disease susceptibility, general health and survival (Strange et al., 2016;Francella et al., 2022;Ratsika et al., 2023).
Research also shows the importance of glucocorticoids (GC) in gut-brain communication, with increased GC levels affecting microbial diversity and composition (Petrullo et al., 2022).Various factors, including, predation (Mohring et al., 2023), heat stress (Lin et al., 2006), high energetic demands during reproduction (Fletcher et al., 2015), and food availability (Romero and Wikelski, 2001) can affect host homeostasis and are known to induce elevated stress levels in animals.Corticosterone, a key glucocorticoid in various vertebrate species, is recognized as a classic endocrine response to stress but also for its role in energy regulation (Almasi et al., 2009).GCs mediate ongoing stress responses, either via maintaining basal levels, allowing other aspects of the stress response to act efficiently, or by actively triggering the stress response (Sapolsky et al., 2000).An alternative view suggests that GCs may suppress the stress response, preventing detrimental over-activation (Sapolsky et al., 2000).
Substantial evidence emphasizes the importance of the microbiota-immune-brain axis in humans and laboratory animals (Fung, 2020).Yet, the understanding of these links in wild animals remains limited (Hird, 2017;Davidson et al., 2020).Studies in birds associated exploratory behavior with microbiota diversity, while learning and memory performance have been correlated with compositional differences (Florkowski and Yorzinski, 2023;Slevin et al., 2020).The gut microbiota's links to stress have also been explored: in common toad tadpoles (Bufo bufo), elevated baseline corticosterone associates with higher microbial diversity (Gabor et al., 2022), while American red squirrels (Tamiasciurus hudsonicus) show lower alpha diversity and fewer gastrointestinal pathogens in response to elevated glucocorticoids (Petrullo et al., 2022).Studies of the microbiota-immune axis in wild barn swallows (Hirundo rustica) and Egyptian fruit bats (Rousettus aegyptiacus) demonstrate that antigen challenges (phytohaemagglutinin and lipopolysaccharides, respectively) can induce changes in gut microbiota composition, which in turn predict the strength of the immune response (Kreisinger et al., 2018;Berman et al., 2023).However, since the majority of these studies have involved adult individuals, it is uncertain whether similar patterns would be observed in nestlings and juveniles.In a rare exception, Stoffel et al. (2020) found that a small proportion of the variation in beta diversity among northern elephant seal juveniles was explained by health status (assessed by counting various white blood cell populations) yet a clear pattern emerged where healthier individuals exhibited higher microbiota diversity.
Here, we investigated the gut microbiota-immune-brain axis in a wild vertebrate population, drawing on predictions based on human and murine models.Using a wild population of common buzzards (Buteo buteo), we sampled nestlings and collected information on bacterial and eukaryotic microbiota diversity (microbiota component).Long-term stress was measured from feather corticosterone (brain component).In addition we performed a series of immune assays (immune component) and estimated the body condition of each individual.Furthermore, we assessed each component at two distinct time points throughout nestling development to incorporate developmental trajectories and dynamic changes.Our approach diverges from traditional methodologies, which often rely on single-point estimates.All components were subsequently integrated into a structural equation modeling framework (Grace, 2006;Lefcheck, 2016) incorporating the following assumptions: 1.All components of the axis exert an influence on the body condition of the nestlings; 2.
The immune system is influenced by both stress levels and the microbiota; 3. Corticosterone levels affect microbiota diversity.In common buzzards we have previously shown that as individuals mature, microbiota diversity declines suggesting a shift toward a stable community following uncontrolled colonization after hatching (Pereira et al., 2024).Building on that knowledge we expect that, as individuals age their gut microbiota diversity will decline, their immune system will mature (increased immune capacity), and stress levels will decrease, leading to a state of homeostasis and improved body condition.

Study system
The common buzzard (Buteo buteo) is a long-lived raptor (up to 25 years) widely distributed across Europe (Walls and Kenward, 2020).This study examines a population that has been closely monitored for over 25 years within a 300 km 2 area (Krüger, 2004;Jonker et al., 2014).Predominantly a resident species, common buzzards breed in tall trees (over 10 meters) from March to July.Adults typically form long-term pairs, build stick nests, and rear on average two nestlings per brood, with reproductive success closely tied to vole availability, their primary food (Kostrzewa and Kostrzewa, 1990;Panek, 2021).During their ≈5-week nestling period (≈35 days), nestlings are entirely dependent on parental care and feeding.Even after fledging (≈47 days), they stay close to the nest and continue receiving parental support (Fachi et al., 2019).Parental care (mostly by the female) decreases during this period: initially, females provide active brooding beyond just feeding times (0-8 days), but as the nestlings grow, female presence decreases (9-30 days), eventually mirroring the care level of males (Hubert et al., 1995).While diet seems to remain consistent (as inferred from observational data), there is a shift in feeding practices: parents transition from directly provisioning nestlings to depositing prey in the nest, allowing the nestlings to gradually consume whole prey items.
We used existing microbiome data from Pereira et al. (2024) and complemented these with new data on feather corticosterone and a set of immune markers.We successfully obtained all sample types from a set of 43 common buzzard individuals, each sampled at two different time points during their nestling phase, totalling 86 samples.Individuals were sampled from 23 nests with an average brood size per nest of two nestlings (Table 1; see Supplementary Tables S1, S2 for n • of ASVs per marker).These were sampled across two different habitats: North of the Teutoburg Forest and south of the Teutoburg Forest (8 • 25 ′ E and 52 • 6 ′ N; Eastern Westphalia, Germany), as described in Krüger (2004) and Jonker et al. (2014).

. Sample collection
In brief, body weight and wing length were recorded at each sampling event, with a nine-day interval on average.The first sampling occurred on average at 19 days of age (mean ± s.d.= 19.3 ± 5.29 days) and the second at 28 days (mean ± s.d.= 27.8 ± 5.16 days).Due to difficulties in precise age estimation prior to sampling (the hatching date was estimated by the number of droppings on the ground below the nest; age was determined following the initial visit to the nest and subsequent measurement of wing length), it was impossible to sample individuals at exactly the same age at sampling points 1 and 2. Nestling age was calculated (postfirst sampling) using a sex-specific polynomial regression model on wing length (Bijlmsa, 1999).Body condition index (BCI) was determined by extracting the residuals of a logarithmic regression of weight on wing length, adjusting for sex.Blood samples (500 µl) were collected from the ulnar vein and stored in heparinized tubes.A small blood drop was used for smears, a portion was stored in ethanol for sex determination, and the remaining volume was centrifuged in the field until separation of plasma was visible.Separated plasma was transferred to a new tube, stored in dry ice, and subsequently stored at −80 • C until further analysis.The remaining red blood cells were resuspended in PBS solution and stored at −20 • C. Cloacal swabs for gut microbiota analysis were obtained and stored in RNAlater, first in dry ice and then long term at −80 • C. In order to assess corticosterone levels, one interscapular feather was pulled from each bird and these were individually stored in paper envelopes.
. Microbiome profiling: DNA isolation, sequencing, and data processing For detailed procedures on microbiome sequencing data processing see Pereira et al. (2024).
Contaminants were identified and removed with the decontam package version 1.18 (Davis et al., 2018).ASVs assigned to Mitochondria, Chloroplast, Vertebrata, Eukaryota, and unassigned taxa were filtered out using QIIME2.Singletons and samples with a minimum frequency below 500 reads were removed.The remaining ASVs were aligned with MAFFT (Katoh et al., 2002) and used to construct a phylogeny with FastTree (Price et al., 2010), both implemented in QIIME2.

. . S rRNA gene sequence data processing
Demultiplexed Illumina sequence data were imported into R version 4.2.2 (R Core Team, 2022).Locus-specific primers were removed using Cutadapt (version 4.4) (Martin et al., 2018).In QIIME2, quality assessment was performed through visualization of quality plots.ASV inference was conducted using DADA2 following the methodology outlined by Callahan et al. (2016).Sequences were trimmed to eliminate low quality regions and paired-end reads were concatenated.In QIIME2, taxonomy was assigned using the naive Bayes taxonomic classifier trained on the SILVA 138.1 database.The decontam pipeline in R was applied to remove contaminants and to perform taxonomy-based filtering (host reads, Mitochondria, Chloroplast, Vertebra and unassigned reads were removed), removal of unique features, and filtering of samples with fewer than 500 reads in QIIME2.The resulting ASVs were aligned using MAFFT, and a phylogeny was constructed using FastTree.

. Assessment of immunity
Due to the complexity of the immune system, we measured four innate immune parameters [bacterial killing activity against Escherichia coli (BKA), lysozyme and haptoglobin concentrations, natural antibodies (HA) and complement (HL) titers] and one component of the acquired immune system [total immunoglobulin Y (IgY) concentration].All these assays are regularly used in wild bird species, including raptor nestlings, both in comparative and within-species studies (see references below).

. . Bacterial killing activity
The bacterial killing activity (BKA) against Escherichia coli (ATCC No 8739) was used to characterize the functional activity of a bird's constitutive innate immune system (Irene Tieleman et al., 2005) using a spectrophotometric version of the assay (Nebel et al., 2021;Brust et al., 2022;Vincze et al., 2022).This assay measures the plasma's capacity to kill microbes ex vivo, determining the organism's ability to eliminate bacterial pathogens encountered, providing an environmentally relevant immune response (Millet et al., 2007;Tieleman, 2018).The assay evaluates the synergic function of several immune components (humoral ones in case of working with plasma), including antibacterial enzymes, complement components, and natural antibodies (French and Neuman-Lee, 2012).Briefly, 12 µl of 1:7 PBS-diluted sample was pipetted in duplicate into 96-well-plate and mixed with 4 µl of ≈1.5 × 105 colony-forming units (CFU)/ml.A positive (not containing any plasma) and a negative control (not containing any E. coli or plasma) was run on each plate.After incubation for 30 min at 37 • C, 83µl of tryptic soy broth (#22092, Fluka) was added to each well.Absorbance at 300 nm was measured with a spectrophotometer (Biotek; µQuant Microplate Spectrophotometer) to determine background absorbance and again after the plates had been incubated for 12 h at 37 • C. The BKA was quantified as the bacteria growth in plasma after 12 h (in %) subtracted by the background absorption in relation to the positive control (Brust et al., 2022).

. . Lysozyme
Lysozyme is an antibacterial enzyme that causes rapid cell lysis, especially in Gram-positive bacteria.It is part of the constitutive innate immune system and is often measured to assess inflammation-induced levels in plasma (Millet et al., 2007).To measure its concentration in plasma, we used the lysoplate assay (Prüter et al., 2020;Brust et al., 2022): 10 µl of sample was inoculated in the test holes of a 1% Noble agar gel (A5431, Sigma) containing 50 mg/100 ml lyophilized Micrococcus lysodeikticus (M3770, Sigma), a bacterium which is particularly sensitive to lysozyme concentration.Crystalline hen egg white lysozyme (L6876, Sigma) (concentration: 0.5, 0.8, 1, 2, 4, 8, 10, 20, and 40 µg/ml) was used to prepare a standard curve for each plate.After 20h incubation at 37 • C, a clear zone developed in the area of the gel surrounding the sample inoculation site corresponding to the bacterial lysis.The diameters of the cleared zones are proportional to the log of the lysozyme concentration.This area was measured three times digitally using the software ImageJ (version 1.48, ImageJ) and the mean was converted to a semi-logarithmic plot into hen egg lysozyme equivalents (HEL equivalents, expressed in µg/mL) according to the standard curve (Prüter et al., 2020;Brust et al., 2022).

. . Haptoglobin
Haptoglobin, is an acute-phase protein in birds and is part of the induced innate immune system.Acute phase proteins are key indicators of immunological function.Their concentrations fluctuate over time, reflecting changes in health and physiological condition (Hõrak et al., 2002(Hõrak et al., , 2003)).Levels can rise quickly in response to infection, inflammation, or trauma (Millet et al., 2007).Elevated plasma haptoglobin often signifies the onset of a non-specific immune response (Matson et al., 2012).We measured haptoglobin concentrations with a commercial kit (TP801, Tri-Delta Diagnostics, Inc.) following the instructions of the manufacturer.Haptoglobin concentrations (mg/ml) in undiluted plasma samples were calculated according to the standard curve on each plate (Prüter et al., 2020;Nebel et al., 2021).

. . Haemolysis-haemagglutination assay
The levels of the natural antibodies and complement were assessed by using a haemolysis-haemagglutination assay as described by Matson et al. (2005), Nebel et al. (2021), andBrust et al. (2022).Natural antibodies, quantified as Haemagglutination (HA) titers, bind non-specifically to various antigens and play crucial roles in opsonization (Ochsenbein et al., 1999;Forthal, 2014).Haemolysis (HL), facilitated by the complement system, is a component of the innate immune system, its activation leads to cell lysis, particularly of cells that have been opsonized (Nauta et al., 2004).After pipetting 10 µl of plasma into the first two columns of a U-shaped 96-well microtiter plate, 10 µl sterile PBS was added to columns 2-12.The content of the second column wells was serially diluted (1:2) until the 11 th column, resulting in a dilution series for each sample from 1/1 to 1/,1024.The last column of the plate was used as negative control, containing PBS only.Ten microliter of 1% rabbit red blood cells (supplied by Innovative Research) suspension was added to all wells and incubated at 37 • C for 90 min.After incubation, in order to increase the visualization of agglutination, the plates were tilted at a 45 • angle at room temperature.Agglutination and lysis, which reflect the activity of the natural antibodies and the interaction between these antibodies and complement (Matson et al., 2005;Prüter et al., 2020), was recorded after 20 and 90 min, respectively.Haemagglutination is characterized by the appearance of clumped red blood cells as a result of antibodies binding multiple antigens, while during haemolysis, the red blood cells are destroyed.Titers of the natural antibodies and complement were given as the log2 of the reciprocal of the highest dilution of serum showing positive haemagglutination or lysis, respectively (Matson et al., 2005;Prüter et al., 2020;Brust et al., 2022).

. . Total immunoglobulin Y concentration
Total IgY, the avian equivalent to mammalian immunoglobulin G, is the primary humoral effector of the adaptive immune system, playing a critical role in neutralizing pathogens (Warr et al., 1995).IgY is essential for long-term immunity, providing protection against recurrent infections.Total IgY was measured using an ELISA with commercial anti-chicken antibodies (Bourgeon et al., 2010;Hanssen et al., 2013).Ninety-six well high-binding ELISA plates (82.1581.200, Sarstedt) were coated with 100 µl of diluted plasma sample (1:4,000 diluted in carbonate-bicarbonate buffer, in duplicates) and incubated first for 1 h at 37 • C and then overnight at 4 • C.After incubation, the plates were washed with a 200 µl solution of PBS and PBS-Tween, before 100 µl of a solution of 1% gelatine in PBS-Tween was added.Plates were then incubated at 37 • C for 1 h, washed with PBS-Tween and 100 µl of polyclonal rabbit antichicken IgY conjugated with peroxidase (A-9046, Sigma) at 1:250 (v/v) was added.Following 2 h incubation at 37 • C, the plates were washed again with PBS-Tween three times.After washing, 100 µl of revealing solution [peroxide diluted 1:1,000 in ABTS (2,20-azinobis-(3-ethylbenzthiazoline-6-sulphonic acid))] was added, and the plates were incubated for 1 h at 37 • C. The final absorbance was measured at 405 nm using a photometric microplate reader (µQuant Microplate Spectrophotometer, Biotek) and subsequently defined as total serum IgY levels (Brust et al., 2022). .

Feather-corticosterone determination
The entire feather was weighed and placed in a tube.For every 10 mg of feather, 1 ml of 100% methanol (p.a.) was added, followed by crushing with scissors.Subsequently, the samples were subjected to an ultrasonic bath incubation at 30 • C for 30 min, followed by overnight incubation in an overhead shaker (25 rpm).On the following day, centrifugation (2 min, 19,800 × g) was performed, and the supernatant transferred to a new tube.The feather pellet underwent two washes, each with half of the extraction volume (2 min each at 19,800 × g); the resulting supernatants were then combined and one further centrifugation was performed for 10 min at 19,800 × g.The supernatant was then filtered through a PTFE membrane (Merck Millipore: Ultrafree-CL), and the filtrate was temporarily stored at −20 • C. A defined amount of methanol (500 µl) from the temporarily stored filtrate was evaporated in a vacuum concentrator.The samples were re-suspended in 250 µl ELISA buffer (Cayman Chemicals Inc. [# 501320]) and stored at −20 • C for at least overnight before conducting the ELISA.Corticosterone concentrations were determined in triplicate following the manufacturer's instructions.f-CORT assay validation is presented in Supplementary Data Sheet 6.

. Statistical analysis
In order to evaluate sequencing depth and sample coverage, rarefaction curves were constructed in QIIME2.Rarefaction was then applied to the 16S rRNA dataset at 4,000 reads and to the 28S rRNA dataset at 2,000 reads.Utilizing the q2-diversity alpha plugin, three alpha diversity metrics were calculated: Shannon diversity index (Shannon, 1948); Faith's Phylogenetic Diversity (Faith PD) (Faith and Baker, 2006) and number of observed ASVs.Prior to constructing the Structural Equation Model (SEM), we explored associations of the variables of interest (BCI, f-CORT, immune assays) with host and environmental factors examined by Pereira et al. (2024).These included sex, habitat and rank (dominance hierarchy within the brood).We have previously established that none of these variables affects microbiota alpha diversity (Pereira et al., 2024).This preliminary analysis aimed to assess the relevance of these variables for the construction of the SEM.For each variable of interest, individual linear mixed models (LMMs) with a Gaussian distribution were constructed using the lmer function from the lme4 package in R (Bates et al., 2015).The significance of factors was assessed through analysis of variance (ANOVA).To account for repeated samples and nest sharing, Nest ID and Individual ID were incorporated as nested random effects (Individual ID nested within Nest ID).Additionally, the Benjamini-Hochberg method was applied to correct the p-values for multiple hypothesis testing (Benjamini and Hochberg, 1995).Scripts and results are presented in Supplementary Data Sheet 2.

. . Structural equation modeling approach
In order to tackle the complex interactions within the gut microbiota-immune-brain axis and their impact on individual condition, we employed a SEM strategy.SEM is a statistical methodology that allows for the simultaneous testing of complex relationships among multiple variables.Integrating factor analysis and multivariate regression analysis, SEM provides enhanced flexibility, enabling variables to both depend on and influence other variables.Furthermore, SEM allows for the incorporation of mediation effects, facilitating the quantification of direct, indirect, and total effects (Grace, 2006;Lefcheck, 2016).Of particular relevance to this study is SEM's capability to construct and model latent variables: ones that cannot be directly measured but are hypothesized to exist (Bollen, 2002;Lefcheck, 2016).

. . Exploratory factor analysis
We started by constructing a latent variable representing Immunity; for this, an exploratory factor analysis (EFA) was performed using the "lavaan" (version 0.6-16) package (Rosseel, 2012) in R. The latent variable model incorporated all immune assays as indicators, with "cluster= Individual ID" set to address repeated sampling and Full Information Maximum Likelihood (FIML) utilized for handling missing data.Estimation was performed using Maximum Likelihood with robust standard errors (MLR), and standardization of values for the latent variable was implemented ("std.lv=T").Model fit was evaluated through the examination of fit indices, including the chi-squared pvalue, Comparative Fit Index (CFI), Root Mean Square Error of Approximation (RMSEA), and Standardized Root Mean Square Residual (SRMR).The conventional "rule of thumb": CFI > 0.9; RMSEA < 0.08; SRMR < 0.08 was used for model fit evaluation (Hu and Bentler, 1999;MacCallum et al., 1996;Sharma et al., 2005).Predicted values for the latent construct were extracted using the "lav.predict"function from "lavaan" with default settings.During the initial assessment, it became evident that haptoglobin concentration did not significantly load on the latent variable Immunity (Figure 1), leading to lower model fit (Table 2).As a result, haptoglobin concentration was subsequently removed and treated individually.

. . Bayesian SEM
A Bayesian structural equation model (Kaplan and Depaoli, 2012) was fitted using the brms R package (Bürkner, 2017;Bürkner, 2018), which incorporates the probabilistic coding language Stan (Carpenter et al., 2017).Four distinct models were inputted into the multivariate analysis in order to test our main hypotheses: Path 1:

BCI ∼ alpha diversity + f-CORT + Immunity + Age + (1|Nest ID/Individual ID)
Path 2: Path 3: Path 4: Path 1 allowed us to infer the contributions of the different components of the axis for nestling body condition, while the remaining paths were designed to capture the dynamics of the gut microbiota-immune-brain. Statistical constraints impose that we assume directionality a priori, although we are aware that these are bidirectional relationships.To account for repeated samples and nest sharing, a nested random effect (individual ID nested within nest ID) was incorporated into the model.Given our sampling design, where individuals of different ages were sampled at each time point (age as a continuous variable), the only way to account for age effects was to include age as a fixed effect.The model was run with four chains, each run with 100k iterations, a warm-up phase of 50k iterations and default priors.Model fit was assessed by examining the convergence of the runs, mixing of chains and performing posterior predictive checks (comparing predicted vs. observed posterior distribution).Marginal and conditional Bayes R 2 (Gelman et al., 2019) were calculated using the "bayes R2" function in brms.Scripts and intermediate results can be found in Supplementary Data Sheets 3, 4.

. Di erential abundance analysis
Analysis of Compositions of Microbiomes with Bias Correction 2 (ANCOM-BC2) with default parameters was implemented in the R package ANCOMBC version 2.0.2 (Lin and Peddada, 2020;Lin et al., 2022).Body condition index, f-CORT, Immunity and age were specified as fixed effects.A nested random effect to account for Individual ID within Nest ID was fitted with the option "rand formula".The Holm-Bonferroni method with a significance cutoff of padj < 0.05 was used to correct P-values for multiple testing (Holm, 1979).Detailed scripts and results are presented in Supplementary Data Sheet 5.

Results
We examined the potential contribution of sex, habitat and brood rank to the variation of BCI, f-CORT, and the various immune assays.The analysis revealed no substantial evidence linking body condition, corticosterone levels, or immune capacity to these variables.However, slight variations in lysozyme and haptoglobin levels were observed among different habitats.These patterns appear to be largely influenced by imbalances in the sample design and the presence of outliers (Supplementary Figure S9).Thus, we assume that the tested variables should not play a significant role in the SEM construct.For a comprehensive overview of the analysis pipeline and detailed results, see Supplementary Data Sheet 2.  .

Exploratory factor analysis
Exploratory factor analysis (EFA) derived a latent variable representing the underlying structure of immune parameter values.We synthesized the results of the immune assays into a composite variable denoted as Immunity, which included all of the immune parameters except haptoglobin concentration (see above).Components of the innate immune system (Hemolysis and Hemagglutination) displayed the strongest factor loadings (λ HA = 0.94, λ HL = 0.90), followed by IgY (λ IgY = 0.50) representing the adaptive immune system (Figure 1).correlated with the studied variables (Supplementary Figure S5).In contrast, when incorporating the Haptoglobin immune assay into the models, no deferentially abundant taxa were identified (Supplementary Figure S6; Supplementary Table S22).

Discussion
The gut microbiota-immune-brain axis shapes a variety of physiological responses through multidirectional communication among the gut microbiota, immune system, and the central nervous system (CNS) (Sylvia and Demas, 2018).This feedback system is known to influence immune function, neuroendocrine pathways, and behavior (Martin et al., 2018;Cryan et al., 2019).Once established, disruptions to this axis may lead to the development of disorders (e.g., irritable bowel syndrome; ulcerative colitis; Alzheimer's and Parkinson's disease in humans) significantly impacting individual health (Rhee et al., 2009;Ghaisas et al., 2016).Here, we used gut microbial diversity measurements, various innate and adaptive immune markers and feather-corticosterone levels (a proxy for stress) to explore relationships among the three components of the axis and the repercussions to body condition in a wild population of common buzzard nestlings.

. Gut microbial diversity, not associated with Immunity and stress
In contrast to our assumptions, our results revealed no association between microbial diversity, Immunity and stress (f-CORT) (Figures 2, 3).However, several parts of the feedback cascade appeared to be functional (see below).Despite the wellestablished influence of the gut microbiota on immune system modulation and its impact on the CNS, these tripartite connections remain elusive in wild populations (Hird, 2017;Davidson et al., 2020;Madden et al., 2022;Pereira et al., 2023), particularly during the early stages of development, in contrast to the studies in humans and mice (Martin et al., 2018;Francella et al., 2022;Lynch et al., 2024).
Contrary to our findings, recent research by Francella et al. (2022) demonstrated clear links between gut microbiota, the immune system and stress in the early-life stages of laboratoryreared mice.Their study showed that immunocompromized mice had increased stress levels, decreased microbial diversity, and alterations in gut microbiome composition post-weaning; additionally, stress impacted the abundance of specific taxa that, in turn, were associated with specific behavioral traits.Furthermore, these behavioral changes were observed after just 1 week of age, demonstrating that stress can interact with host immunity  during early development.Discrepancies between the results of the mice study and ours can be explained by various non-mutually exclusive factors.It is important to consider that their research was conducted under laboratory conditions, which lack ecologically relevant environmental factors, whereas our study was performed under ecologically realistic, natural circumstances.Indeed captivity is known to have an effect on each components of the axis (Slevin et al., 2020;Florkowski and Yorzinski, 2023), the studies also focused on different aspects of the immunity (cellular vs. humoral) and made use of immunocompromised/gene knockout mice.Furthermore species have different diets and life-history strategies in terms of fast-slow life-history continuum (Réale et al., 2010;Bing et al., 2022).
Despite the scarcity of such distinct gut microbiota-immunebrain studies in other vertebrate species with lower degrees of experimental manipulation, studies on captive organisms have demonstrated a connection between the brain and the gut microbiota.For instance, higher gut microbiota alpha diversity in zebra finches (Taeniopygia guttata) has been correlated with elevated exploratory behavior, and in house sparrows (Passer domesticus), beta diversity was associated with enhanced cognitive performance (Slevin et al., 2020;Florkowski and Yorzinski, 2023).It has been recognized that the gut microbiome modulates stress responses, particularly through the hypothalamicpituitary-adrenal (HPA) axis (Foster and Neufeld, 2013); in general, higher stress levels tend to be associated with reduced microbiota diversity, specifically a reduction in the amount of rare and pathogenic taxa (Petrullo et al., 2022).While the notion of diverse microbiota being crucial for robust immunity is widely accepted (Hooper et al., 2012;Lozupone et al., 2012), extending these conclusions beyond a small number of experimental study systems is not straightforward, as several studies have found no association between microbiota diversity and immune indices.Instead, they observed that compositional differences or changes in specific taxa were more closely linked with variations in immune markers (Kreisinger et al., 2018;van Veelen et al., 2020;Fleischer et al., 2024).Our focus on microbial diversity, without considering composition (beta diversity), limits the depth of our analysis (Shade, 2017;Reese and Dunn, 2018).However, incorporating the complexities of multivariate metrics into structural equation modeling remains challenging.Nevertheless, we hope that ongoing advancements in methodologies for analyzing compositional data (Sweeny et al., 2023;Fountain-Jones et al., 2024) will facilitate this approach in future studies.It is also important to consider timing in the context of stress responses, as short-term and prolonged stress can induce contrasting physiological reactions (Martin, 2009).Evaluating different measures of HPA axis activity would offer a more comprehensive understanding of stress response dynamics (Stothart et al., 2019).Short-term acute stress triggers an organism's immune response (Martin, 2009) and enhances intestinal mucus secretion (Castagliuolo et al., 1996).Conversely, prolonged stress typically suppresses immunity and reduces mucus production, impacting the microbiota differently (Estienne et al., 2010;Shigeshiro et al., 2012).
Several studies of wild populations have delved into some components of the gut microbiota-immune-brain axis (Noguera et al., 2018;Petrullo et al., 2022;Berman et al., 2023), yet few have comprehensively addressed the entire system.An exception is a recent study on eastern newts (Notophthalmus viridescens), which explored the effects of experimentally elevated CORT levels on various immune indices and skin microbiome but did not find any clear evidence for a relationship among the three components (Pereira et al., 2023).
The examples provided demonstrate the challenges encountered in establishing associations among the components of the gut microbiome-immune-brain axis, either due to the complexities of wild settings or the specific time window investigated, or indeed, because there are no discernible associations in the corresponding systems.These challenges increase in less controlled study systems (further aggravated by the fact that studies of completely wild populations are almost non-existent), particularly when dealing with systems where certain components of the axis lack experimental challenge, alteration, or complete knockout.Given this, and considering the absence of an evident connection with the microbiome axis, we propose the following possible, not mutually exclusive mechanisms to be at play: 1.The gut microbiota may not exert a significant influence in the initial stages of life, especially when compared to the impacts of stress and the immune system.This is supported by the robust link between f-CORT and Immunity, which we find despite the correlative nature of our study system (Figures 2-4A).Consequently, we suggest that stress and immune regulation may hold greater importance for maintaining homeostasis during early development.2. The complexity of wild environments introduces many background effects, potentially including diverse diets, exposure to various pathogens, and the unpredictable nature of ecological interactions.These factors might overshadow the subtle and context-dependent relationships observed in more controlled settings.In reality, the connections found in these controlled environments may not be as crucial or representative of the complete picture encountered in the more realistic complexities of a wild setting.3. Direct comparisons between species and generalizations across species might be challenging due to ecological differences between the systems.

. Stress, Immunity and body condition
The relationship between stress, immunity, and body condition is a strong feedback mechanism in animal physiology (Vagasi et al., 2018).Baseline levels of corticosterone play important roles in metabolism, development, reproduction, behavior, and immunity (Sapolsky et al., 2000).While beneficial in the short term for resolving inflammation and preventing an overshoot of the immune responses (Dhabhar, 2018), prolonged stress-induced elevation of corticosterone compromises the immune system over time (Dhabhar and Mcewen, 1997).A study of house sparrows revealed that prolonged activation of the stress response inhibits components of the innate immune system, such as complementmediated lysis, bacteria-killing ability, and agglutination.Indeed, numerous studies have demonstrated the detrimental effects of high stress levels, for prolonged periods of time, on overall body condition, survival, and reproductive success (Angelier et al., 2010;Mikkelsen et al., 2023;Quirici et al., 2021).
We used feather CORT as an indicator of long-term stress.In birds, circulating CORT accumulates in developing feathers, serving as a cumulative gauge of hormone concentrations during feather growth (Jenni-Eiermann et al., 2015;Romero and Fairhurst, 2016).Our results show a negative association between Immunity and f-CORT, supporting the established concept that immunosuppression is expected in the face of allostatic overload (chronic stress) (Romero et al., 2009).Additionally, prolonged stress can incur costs for developing individuals, may divert resources away from essential physiological processes, and impair an individual's ability to fend off infections and maintain overall wellbeing, as evidenced here by a decline in body condition (Figures 2-4A, B).

. Bacterial microbiota diversity and body condition
Our results show that nestlings with higher Shannon diversity have better body condition (Figures 2, 4C).This positive association can be attributed to the enhanced resistance of more diverse gut communities against pathogenic invasion, increased stability and resilience to disturbance (Buffie and Pamer, 2013).Additionally, a more diverse microbiota can potentially offer a broader range of functions performed by various bacterial taxa, leading to benefits for the host (Heiman and Greenway, 2016).Conversely, lower microbiota diversity is typically viewed as being detrimental to hosts (Le Chatelier et al., 2013), as it implies a loss of essential functions that could result in reduced nutrient assimilation or immunodeficiency (Round and Mazmanian, 2009;Hanning and Diaz-Sanchez, 2015).However, it has also been shown that high diversity might be associated with a state of dysbiosis, so a reduction in diversity could signal a return to homeostasis (Johnson and Burnet, 2016;Kohl et al., 2018;Coyte and Rakoff-Nahoum, 2019).Our results align with the Anna Karenina principle, which states that changes in microbiota due to disturbances lead to shifts from stable to unstable community states (Zaneveld et al., 2017).We propose that the increase in Shannon diversity in buzzard nestlings indicates an increase in stable and abundant taxa, which may offer greater benefits to the host (Hanning and Diaz-Sanchez, 2015).Metrics like Faith's PD and the number of ASVs treat all taxa equally (Faith and Baker, 2006), potentially explaining why they do not show associations with body condition (Figure 3; Supplementary Figure S1).
The absence of clear links between eukaryotic microbiota diversity and body condition may indicate a slower rate of change in the eukaryotic microbiota.Bacterial and eukaryotic taxa are likely to have distinct roles within the gut ecosystem (Oever and Netea, 2014;Vemuri et al., 2020).Fast-changing bacterial taxa may experience stronger competition, leading to positive selection on core functionally-relevant taxa (Abt and Pamer, 2014;Coyte and Rakoff-Nahoum, 2019), whereas changes in the more complex and less abundant eukaryotic taxa (Chin et al., 2020) may occur at a slower pace and therefore not be observed during early stages of host development.
. Age e ects on stress, Immunity and eukaryotic microbiota diversity The positive relationship between age and Immunity, and the negative association of age and stress levels (f-CORT; Figures 2,  3), most likely reflect nestling development and the maturation of the immune system.Initially, nestlings rely on innate and maternal transferred immunity (Klasing and Leshchinsky, 1999;Palacios et al., 2009).As they mature, exposure to pathogens and foreign microbes through the environment and diet, contact with nest material and siblings, challenges and stimulates the immune system, thus aiding its maturation (Lochmiller and Deerenberg, 2000;Morais et al., 2020;Oldereid et al., 2023).During the ontogenetic development of birds, exposure to corticosterone occurs in both the prenatal (embryonic) and postnatal (nestling and fledgling) periods (Henriksen et al., 2011;Strange et al., 2016).In the prenatal stage, corticosterone is transferred from the mother to the embryo through the egg yolk, and is influenced by the maternal environment, e.g., exposure to predators, competitors and other stress-inducing cues during egg production, and differences in environmental quality (Hayward and Wingfield, 2004;Saino et al., 2005;Love et al., 2008).In the postnatal stage, an initial surge in corticosterone levels may serve as an adaptive response to stressors associated with hatching, exposure to a new environment, and nutritional demands, being beneficial in the short term (Chin et al., 2009;Spencer et al., 2009;Strange et al., 2016).As nestlings mature, the HPA axis undergoes maturation, leading to improved stress response regulation.Initially, nestlings invest significant energy in growth and development.As they progress through early life stages, maturation enables more efficient allocation and prioritization of resources (Smulders, 2021;Spencer et al., 2009).As mentioned earlier, the maturation of both physiological systems is not independent; instead, a bidirectional interaction regulates both the immune system and the HPA axis (Francella et al., 2022).
We find a decline in eukaryotic microbiota phylogenetic diversity and the number of ASVs with age.Dominant taxa show higher resilience to disturbance, securing their positions via selective filtering and the occupation of core niches (Costello et al., 2012;Abt and Pamer, 2014;Coyte and Rakoff-Nahoum, 2019).Rare taxa are primarily acquired through stochastic processes but still significantly contribute to certain community alpha diversity measures (like Faith PD and n • of observed ASVs).Over time, these rare taxa are progressively replaced by the more dominant taxa (Shade et al., 2014).It's important to emphasize that microbiota colonization can be either deterministic or stochastic.Studies in gnotobiotic animals highlight this variability, demonstrating that some bacterial strains colonize in a deterministic manner, while others do so stochastically (Vega and Gore, 2017;Jones et al., 2022;Hayashi et al., 2024).Understanding how these processes contribute to community assembly during development in wild populations remains challenging.The dominant taxa also show closer phylogenetic relationships to one another, suggesting potential specialization or competitive advantages driving their prevalence (Janiak et al., 2021;West et al., 2022;Davies et al., 2022).This hints at an initial diverse gut microbiota in newborns due to rapid post-hatching colonization (Trevelline et al., 2018;Pereira et al., 2024).Additionally the decline in microbiota diversity could be linked to decreased parental care (Hubert et al., 1995), resulting in less microbial transmission from parents to nestlings, and to dietary changes as nestlings shift from being fed individual prey pieces by their parents to consuming whole prey items left in the nest (Hoffmann et al., 2013;David et al., 2014).

Conclusion
We investigated the interactions of the gut microbiotaimmune-brain axis in raptor nestlings.As far as we are aware, this represents one of few studies exploring this axis in a wild vertebrate population and incorporating different time points during the nestling phase in a longitudinal study design.
While there was no conclusive evidence for the microbiotaimmune-brain axis in our study system, we did find evidence for the hypothesized relationships among stress, Immunity, and body condition.Elevated f-CORT levels were linked to immunosuppression and adverse effects on overall body condition, suggesting that immune and stress regulation play a dominant role in early nestling development.Shannon bacterial microbiota diversity was positively correlated with nestling body condition, suggesting a potential benefit of a diverse and stable gut microbiota.Age plays a crucial role in influencing immune development, stress responses, and eukaryotic microbial diversity.The decline in eukaryotic microbial diversity with age implies an early uncontrolled gut colonization, followed by selective removal of non-relevant taxa.Our study thereby contributes to a growing body of knowledge on the dynamics of the gut microbiota-immune-brain axis in wild populations.

FIGURE
FIGURELatent variable model for "Immunity".Arrows point to di erent immune assays, which serve as components of the latent variable.Values are indicating the strength of the relationship between the latent variable and its indicators.Red values indicate variables that loaded significantly in the model while black values denote non-significant loadings.*** p < .; * p < . .

FIGURE
FIGUREDiagrams representing the Bayesian structural equation models for the di erent bacterial microbiota diversity metrics.Illustrated are the results of the models incorporating the latent variable, with results for the Haptoglobin assay superimposed.Values represent the estimated mean e ect of a predictor (µ) on the outcome variable.Credible intervals ( % CI) are provided for predictors that explain the outcome variable.

FIGURE
FIGUREDiagrams representing the Bayesian structural equation models for the di erent eukaryotic microbiota diversity metrics.Illustrated are the results of the models incorporating the latent variable, with results for the Haptoglobin assay superimposed.Values represent the estimated mean e ect of a predictor (µ) on the outcome variable.Credible intervals ( % CI) are provided for predictors that explain the outcome variable.

FIGURE
FIGURE Regression plots illustrating the main results of the SEM model.(A) Negative relationship between immunity and f-CORT; (B) Negative correlation between body condition and f-CORT; (C) in bacterial microbiota Shannon diversity with body condition.
TABLE Model fit indices for exploratory factor analysis.
CFI, Comparative Fit Index; RMSEA, Root Mean Square Error of Approximation; SRMR, Standardized Root Mean Square Residual.The presented p-values correspond to the Chi-square test of model fit.